function [ W ] = MyDeconv2D( H,S,I,J, iter)
%MYDECONV2D Summary of this function goes here
%   Detailed explanation goes here
    
    
    [K,L] = size(S);
    W0 = zeros(I,J);
    W1 = W0;
    %W_1 init:
    for i = 1:I
        for j = 1:J
            
            W0(i,j) = GenWij1(i,j,K,L,I,J, H,S);
        end
    end
    for count = 2:iter
        for i = 1:I
            for j = 1:J
                W1(i,j) = W0(i,j)*GenWijN(i,j,K,L,I,J, H,S,W0);
            end
        end
        W0 = W1;   
    end
    W = W0;
end

function y = GenWij1(i,j,K,L,I,J, H,S)
    
    y = 0;
    
    e = i+K-1;
    f = j+L-1;
    
    for m = i:e
        for n = j:f
            x1 = H(m,n)*S(m-i+1,n-j+1);
            
            a = max(1, m-K+1);
            b = min(m,I);
            c = max(1,n-L+1);
            d = min(n,J);
            
            x2 = 0;
            for p = a:b
                for q = c:d
                    x2 = x2+(S(m-p+1, n-q+1));
                end
            end
            
            y = y+(x1./x2);            
        end
    end

end

function y = GenWijN(i,j,K,L,I,J, H,S,W)

    y = 0;
    
    e = i+K-1;
    f = j+L-1;
    
    for m = i:e
        for n = j:f
            x1 = H(m,n)*S(m-i+1,n-j+1);
            
            a = max(1, m-K+1);
            b = min(m,I);
            c = max(1,n-L+1);
            d = min(n,J);
            
            x2 = 0;
            for p = a:b
                for q = c:d
                    x2 = x2+(W(p,q)*S(m-p+1, n-q+1));
                end
            end
            
            y = y+(x1./x2);            
        end
    end

end


